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ABSTRACT 

We perform a Bayesian parameter inference in the context of resonantly damped 
transverse coronal loop oscillations. The forward problem is solved in terms of para- 
_H ■ metric results for kink waves in one-dimensional flux tubes in the thin tube and thin 

Ph' boundary approximations. For the inverse problem, we adopt a Bayesian approach to 

O ■ infer the most probable values of the relevant parameters, for given observed periods 

-^ ■ and damping times, and to extract their confidence levels. The posterior probability 

, "^z distribution functions are obtained by means of Markov Chain Monte Carlo simu- 

lations, incorporating observed uncertainties in a consistent manner. We find well 
^ ■ localized solutions in the posterior probability distribution functions for two of the 

£r} ■ three parameters of interest, namely the Alfven travel time and the transverse inhomo- 

OK ■ geneity length-scale. The obtained estimates for the Alfven travel time are consistent 

. . with previous inversion results, but the method enables us to additionally constrain 

Q ■ the transverse inhomogeneity length-scale and to estimate real error bars for each pa- 

rameter. When observational estimates for the density contrast are used, the method 
enables us to fully constrain the three parameters of interest. These results can serve 
to improve our current estimates of unknown physical parameters in coronal loops and 
to test the assumed theoretical model. 
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1. Introduction 



Magnetohydrodynamic (MHD) seismology stands as one of the few indirect methods for the 
determination of difficult to measure physical parameters in solar coronal structures. It relies on the 
combined use of observed oscillatory properties of MHD waves in solar atmospheric magnetic and 
plasma structu res together with theoret ically obtained wave propertie s. It was first suggested b y 
Uchidal (119701) and Roberts et al.l (119841) . in the coronal context, and by iTandberg-HanssenI (| 1995b . 
in the context of solar prominences. The last years increase in the number and quality of observa- 
tions of wave activity in the solar atmosphere and the refinement of theoretical models have allowed 
the practical implementation of this technique. In the context of coronal loops, coronal seismol- 
ogy has allowed the estimation and/or restriction of releva nt parameters such as the magnetic field 



strength (INakariakoy & OfmanI 120011) . the Alfven speed (IZaaarashvili 



Goossens et al. 



20081). the transversal densitv structuring (IVerwichte et al. 



2003 



Arregui et all 12007 



20061), or the coronal 



density scale height (|Andries et al.ll2005l : IVerth et al.ll2008l) . Of particular r elevance has been th e 
use of the concept of p eriod ratios as a se ismological tool, firs t pointed out by lAndries et al.l(|2005|) : 
Goossens et al.l (|2 



20061) . and reviewed bv lAndries et all (I2009I) . 



Most of the efi'orts in this area have been concentr ated on the phenomenon o f qui ckly damped 
transv erse oscillations in coronal loops, first reported by lAschwanden et al.l (Il999() and lNakariakov et al 
(Il999|) . These oscillations are interpreted in terms of linear MHD kink waves of a magnetic flux 
tube, a wave mode with mixed fa st and Alfven charact er, its Alfvenic nature being dominant in 
and around the resonant position (IGoossens et al.ll2009|) . where the global eigenmode frequency 
matches the local Alfven frequency. Starting with the simplest model that con siders the fun- 
damental transverse oscillation of a magnetic flux tube (lEdwin & Robertslll983h several model 
improvements have included other efl'ects, such as the curvature o f coronal loops, density stra tifi- 
cation, or the departure from circular cross section of the tubes (see lRuderman & Erdelyill2009l for 
a recent review). All these new ingredients have been seen to produce second order eff"ects on the 
main wave properties. 



Th e advancement in seisrnology of kink modes in coronal loops is reviewed by ^ Goossens 



(|2008|) . iNakariakov & OfmanI (1200 If) used the observed periods and theoretical estimates of the 
periods, based on the long wavelength approximation for a uniform coronal loop model, to derive 



estimates for the magnetic field strength, by making assumptions on the density. iGoossens et al. 



(I2OO2I) used the observed damping rates and theoretical values of the damping rates, based on the 
thin boundary approximation, to derive estimates for the rad ial inhomogeneity length- scale, once 
the density contrast was assum ed. lAschwanden et al.l (l2003h vi sed the observed damping rates and 
the damping rates computed by IVan Doorsselaere et al.l (|2004) . outside the thin boundary approx- 
imation, to compare the theoretically predicted density contrasts to estimates obtained indirectly 
taking into account efl^ects such as the subtraction of the background flux, the line-of-sight inte- 
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gration of the emission measure, the spatial smearing due to the transverse motion, and the point- 
spread function of the instrument. In their study, the internal density could only be determined as 
a function of the external density. 

The main limitation of these studies is that the determination of th e magnetic field stren gth 
(INakariakov & Ofmanll200lb or transverse inhomogeneity length-scale (|Goossens et al.l 120021) is 
only possible if a given value for the density contrast is adopted. Otherwise, an infinite number of 
solutions to the i nverse problem arise . Thi s was shown in the nu merical and analytic seismologi- 
cal inversions by Arregu i et al.l (|2007h and lGoossens et al.l (l2008h . that combine the observational 
information on both periods and damping times for resonantly damped kink modes in a consis- 
tent manner. These inversions enable to infer information about both the internal Alfven speed 
and the transverse density structuring. Their approach is to make no assumption on the particular 
value of any of the physical param eters of interest. The ir ranges of variation, compatible with 
observations, are instead obtained (|Goossens et al.ll2008h . Following this approach, a complete 
solution to the inverse problem is obtained, from which general properties and liming cases can 
be studied. The solution gives rise to an infinite number of equally valid equilibrium models that 
explain observations, characterized by three parameters: the density contrast, the transverse inho- 
m ogeneity length-scale , and the Alfven speed. In a ddition, in the inversion techniques presented 
by lArregui et al.l (|2007|) and lOoossens et al.l (|2008[) . it is not straightforward to devise a consis- 
tent method to compute uncertainties on the inferred parameters from the measurement errors on 
observed periods and damping rates. 

Bayesian inversion techniques for parameter inference can help us to overcome these lim- 
itations and are applied in this paper to the determination of physical parameters in oscillating 
coronal l oops. These te chniques have proven to be useful in a number of areas of physical 



scien ces (IGregorv 



2005h. i ncludi ng solar physics. In the context of inversion of Stokes pro- 



files, lAsensio Ramos et al.l (l2007h apply Bayesian techniques to analyze the performance of a 
given radiative transfer model to fit a given observed Stokes vector aiming at obtaining infor- 
mati on about the thermod y namic and rnagneti c properties of solar and stellar atmospheres (see 
also lAsensio Ramosll201 1|) . iMarsh et al.l (|2008|) use a Bayesian probability -based approach to the 
problem of detecting and parameterizing oscillations in the upper solar atmosphere, in order to de- 
termine the nurnber o f oscillations present, and their properties. In the context of coronal heating, 
Adamakis et al.l (|2010) employ Bayesian mod el comparison techn iques to determine the preferred 
heating location along coronal loops. Finally. Ilreland et al.l (120101) have recently presented an au- 
tomated detection technique for oscillating regions in the solar atmosphere based on Bayesian 
spectral analysis of time series and image filtering. 

In this paper, a Bayesian MHD seismology inversion technique is presented and the results 
from its application to coronal loops are described. The layout of the paper is as follows. Section|2] 



describes the physical model and the forward and inverse problems. In Section [31 the details of the 
developed Bayesian technique are given. The analysis and results of its application to synthetic 
data and to real coronal loop oscillations are discussed in Section |4l Finally, in Section [51 our 
conclusions are presented. 



2. The forward and inverse problems 

The application of the Bayesian formalism to parameter inference is rather general and can 
be applied to any model that explains a given set of observations. In our case the theoretical 
model is the resonantly damped MHD kink mode interpretation of quickly damped transverse 
oscillations of coronal loops. The classic theoretical model assumes that coronal loops can be 
represented as straight cylindrically symmetric magnetic flux tubes with a uniform magnetic field 
pointing along the axis of the tube. In the zero plasma-jS approximation coronal loops are density 
enhancements with a constant internal density, pi, a constant external density, pe < Pi, and a 
non-uniform transitional layer of thickness I that connects both regions. Analytical theory for 
linear MHD kink oscillations based on t he thin tube and thin boundary approximations gives us 
two equations (see iGoossens et al.ll2008[) for the period and damping time of resonantly damped 
oscillations. 



P = .„V2f£iir and Ii = i£±ii. 
I ^ / P n^-ll/R 



(1) 



The factor 2/n in the damping rate expression has its origin on the assumption of a sinusoidal 
variation of density across the non-uniform layer. These equations express the period, P, and the 
damping time, r^, which are observable quantities in terms of the internal Alfven travel time, tai, 
the density contrast, ^ = pi/pe, and the transverse inhomogeneity length scale, l/R, in units of the 
radius of the loop. These three quantities {tai, ^, and l/R) are the seismic quantities in the sense 
that they are the quantities that we aim to determine with the use of observed values for the period 
and the damping time. Since we only have two equations that relate the three unknown quantities 
to the two observe d quantities there is an infinite number of solutions to the inverse problem, as 
first pointed out by lArregui et al.l(|2007|) . These solutions have to follow a precise one-dimensional 
solution curve in the three-dimensional parameter space of unknowns. This curve constitutes a 
characteristic a ttribute of resonant a bsorption subject to observationa l testing. It was first ob tained 
numerically by lArregui et al.l(|2007|) and subsequently analytically by iGoossens et al.l(|2008|) under 
the thin tube and thin boundary approximations. 

Because of the ve ry good agreem e nt bet ween the analytic and numerical inversion solution 
curves demonstrated in IGoossens et al.l (l2008b . in the following we adopt the analytic approxima- 



tion for the forward problem. The clear advantage is that the forward problem becomes algebraic 
and this simplifies the statistical inversion scheme considerably. For a particular observed event, 
with fixed period and damping time, the seismic variables are constrained to the following intervals 
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with C = TTTd/P known from observations. The seismic variables are no t indep endent, but are 
related to one another by the following functions defined in lGoossens et al.l (|2008|) 
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Of the four functions only two are independent. This allows to consider one of the three unknowns 
as a parameter, say z = 1/2R, let it take on all values in /^ and then compute the corresponding 
values of y = F](G2(z)) and 4" = ^2(2). An example of the resulting inversion curve is shown in 
Figure 1 in lGoossens et al.1 dlOOSh (see also Figure [2l in this paper). 



Although all the possible combinations of parameters along the seismic curve equally well 
explain observations it turns out that the obtained curve gives rather constrained values for the 
Alfven travel time. Moreover, when applied to prominence thread transverse oscillations, given 
their large density contrast values, the asymptotic behavior of the solution curve towards large 
density contrast values allows to derive asymptotic values for the internal Alfven speed and the 
transverse inhomogeneity length scale (|Arregui & Ballesteiil2010|) . In spite of these virtues, the 
problem with the infinite number of solutions remains. A genuine measurement of any of the 
three unknowns would enable us to collapse the one-dimensional inversion curve into a zero- 
dimensional solution point. The matter is to obtain reliable observational estimates for any of the 
three unknowns. Of related importance is to devise a consistent procedure for incorporating to the 
inversion schemes measurement errors in periods and damping times to compute the uncertainties 
in the inferred physical parameters. 
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3. Probability model and numerical method 

3.1. Mathematical formulation 

Our model describing resonantly damped transverse loop oscillations contains three free pa- 
rameters, ^, l/R, and TAi that have to be estimated. Let these three parameters be gathered in 
the vector 6=[^, l/R, tai]' and let the set <i=[/'obs,Td,obs]' contain discretized data on period and 
damping times. When data, as a result of measurements determine th e values of d, the sta te of 



knowledge on 6 is represented by p(6\d), which is given by Bayes' rule (|Bayes & Pricelll763|) 



p{0\d) = — — . (3) 

P(d) 

In this expression, p(0\d) is the so-called posterior probability distribution. The function p{d\d) is 
the conditional probability distribution of the data given the parameters. It contains the theoretical 
relations between parameters and data (including the noise properties) and is a measure of how 
well the data are predicted by the model. Before d has been observed, it represents the probability 
distribution associated with possible data realizations for a fixed parameter vector. A posteriori, 
after observation, p(d\0) has a very different interpretation. It is the likelihood of obtaining a real- 
ization actually observed as a function of the parameter vector, and is hence called the likelihood 
function. Under the assumption that observations are corrupted with Gaussian noise and that they 
are statistically independent, the likelihood can be expressed as 



p{d\6) = {IncrpCTr) ^ exp 



vp - p'^'^mf ^ [^d-rTiO)^ ^^^ 



lo-l lal 



with P^^^(6) and t]^{6) given by the forward analytical problem of Equation ([T]). Likewise, cr], and 
cr\ are the variances associated to the period and damping times, respectively. 

The quantity p{6) is the prior probability of the model vector. It should contain any prior 
information we might have on the model parameters, without taking into account the observed 
data. Three different versions are used in this paper. The normalization constant on the right- 
hand side of Equation ([3]) is the marginal likelihood, called "Bayesian evidence" in cosmology 
applications. The evidence can be written as 



= / 



p{d) = I pid\e)pie)de, (5) 

which is an integral of the likelihood over the prior distribution that normalizes the likelihood 
and turns it into a probability. This quantity is central for model comparison purposes, since it 
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evaluates the model's performance in the light of data. In Bayesian inference, to compute the 
posterior distribution one needs to run the model over the entire parameter space of interest and 
sum that all up. 

In the Bayesian formulation the inference of a given parameter set 6 is based on the knowledge 
of p(6\d), which contains all the information available on 6 given the data d and therefore is the 
solution to the inverse problem. Once the posterior distribution p{6\d) is known, the position of the 
maximum value gives the most probable combination of parameters 9 that fit the data. Bayes' rule 
constitutes a simple mathematical formulation of the process of learning from experience. What 
can be inferred about the model parameters "a posteriori" is a combination of what is known "a 
priori", independently of the data, and of the information contained in the data. It is an appealing 
formulation that shows how previous knowledge can be updated when new information becomes 
available. An additional advantage of the Bayesian approach is that error propagation is consis- 
tently taken into account by marginalization. When one needs to know how the data d constrains 
one of the parameters, say 6i, it is enough to compute the following integral 



= / 



p{ei\d) = pie\d)dei . . . dOi-ydeM . . . dON. (6) 

The result is known as the marginal posterior, which encodes all information for model parameter 
6i available in the priors and the data. Furthermore, it correctly propagates uncertainties in the rest 
of parameters to the one of interest. 



3.2. Prior information 

The prior choice, p(0), is a fundamental ingredient of Bayesian statistics, in order to obtain 
optimal results. The results obtained from the statistical inference should be independent on prior 
information, provided the prior has a support that is non-zero in regions of the parameter space 
where the likelihood is large. In this case, repeated application of Bayes' rule will lead to a pos- 
terior probab ility distribution that converges to a common objective inference on the hypothesis 



(|Trottall2008l) . Two scientists in the same state of knowledge should assign the same prior, hence 
the posterior should be identical if they observe the same data. In our particular case, dealing with 
the inversion of physical parameters in oscillating coronal loops, expressions © define the inter- 
vals over which the unknown physi c al par ameters are allowed to vary, according to the analytic 



inversion scheme by iGoossens et al.l (|2008|) . 



When all information we have on unknown parameters is restricted to ranges of variation, a 
reasonable choice is to assign the same probability to all values contained within that range. This 
defines a uniform type prior, and we can write 
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where H(x, a,b)is the top-hat function 



(7) 



^ -r— a < X < b, 

b-a ' 



H(x, a, b) 



(8) 







otherwise. 



Another option is to assign a decreasing probability distribution for increasing values of the 
parameter, over a given range. This can be accomplished by considering a Jeffreys' type prior, and 
we can write 



pm = 



0,log 



9Y 



1-1 



(9) 



If some additional information about the unknown parameter is available from observations, 
a Gaussian distribution centered on the measured value can be used, so that 



-1/2 



piOi) = (incrll) exp 



■{0-lus.f 



2^1. 



(10) 



In this expression, the mean //g. and the standard deviation cr^. would be directly obtained from 
observations. 

In our analysis, we use a uniform prior for the transverse inhomogeneity length-scale and 
the internal Alfven travel time. The limits for the transverse inhomogeneity are imposed by the 
modeling, since this parameter must be in the range l/R 6 [0, 2]. For the internal Alfven travel time, 
we have considered uniform prior distributions with minimum and maximum values that enclose 
the intervals given by /y in expressions ^, taking into account the observed period. Coronal loops 
are over-dense structures with an internal density that is at most one order of magnitude larger that 
the external coronal density. Densities two orders of magnitude larger than the coronal density, 
typical in prominence plasmas, are less likely in coronal loops. Our analysis makes use of the 
three prior distributions given by Equations ([7]), dH), and (flOl) . for the density contrast. The virtues 
and disadvantages of each of them are discussed in Section HTTl Figure[I]displays an example of the 
three different types of prior in density contrast. In each case, the integral of the prior probability 
over the entire range of values must amount to one. 



3.3. Markov Chain Monte Carlo method 

The objective of our parameter inference under the Bayesian framework is to sample the full 
posterior distribution so that, for instance, we can locate the most plausible model that maximizes 
Equation © or calculate the marginalization integrals. To this end, one has to evaluate p{0\d) for 
different combinations of parameters. If ten values for parameter are to be evaluated and Npa^ is the 
number of parameters this means that ~ lO^p" evaluations are needed. This exponential increase 
of the number of evaluations with the number of parameters is known as the curse of dimension- 
ality. A convenient way to handle this kind of problem is to perform a Markov Chain Monte 
Carlo (MCMC) simulation, for the sampling of the posterior probability distribution function. The 
MCMC algorithm we apply allows us to construct a sequence of points in parameter space, called 
a "chain", whose density is proportional to the posterior distribution function. It therefore provides 
a method for sam pling the posterior distributio n up to a multiplicative constant. We use the same 
code employed byl Asensio Ramos et al.l (|2007 h in their inversion of Stokes profiles. The details of 



the technique and the algorithm are given in that paper, so only the most relevant information is 
detailed here. 

The obtained Markov chain is defined as a sequence of random variables (^o, G\, . . ., d„), 
such that the probability of the Oi element in the chain only depends on the value of the previous 
element, 6i_i. Starting from a given element in the chain, the next point is chosen in such a way that 
the distribution of the chain asymptotically tends to be equal to the posterior distribution. Markov 
chains can be shown to converge to a stationary state where successive elements of the chain are 
samples from the target distribution, in our case the posterior p{6\d). 



Our method uses the Metropolis algorithm (IMetropolis et al.ll 19531 : lNeallll993h . Starting from 



a given vector of parameters 6q, the posterior probability given the data is calculated, p{6o\d). 
This includes the calculation of the priors and the likelihood, which involves the evaluation of 
the forward problem. Next, a new vector of parameters, Ot, is obtained sampling from a proposal 
density distribution, q{6i\6i-{). Again the posterior is evaluated, p{0i\d). The ratio 

pmd)qmOi-i) 



p(ei_i\d)q(ei_im 

is evaluated and 6i is admitted with probability /3 = min[l, r]. The process further progresses by 
obtaining a new vector of parameters from the proposal density distribution. 

The key ingredient in this process is the proposal density q(6i\6i-{). This distribution is used 
to sample points from the posterior starting from a given point. Ideally, it should be chosen as 
close as possible to the unknown posterior distribution. Also, it should be simple to sample from 
it. Since the first condition is impossible to fulfill unless the posterior is known, we choose a suf- 
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ficiently general distribution function which is easy to sample from. The obvious selection is a 
multivariate Gaussian distribution. In our case, we allow for a non-diagonal covariance matrix in 
the Gaussia n distribution, which impro ves the sampling efficiency of the MCMC code (see Ap- 



pendix A in lAsensio Ramos et al.ll2007l for more details). We have chosen a uniform distribution 
for the initial A^unif steps of the chain. Once some information about the posterior is known the 
algorithm changes to a Gaussian proposal density centered on the current value of the parameter. 



4. Statistical inversion results 

4.1. Synthetic data 

In order to assess the performance of our MCMC code, when applied to the problem at hand, 
the inversion of physical parameters is first performed with synthetic data for period and damping 
times. These synthetic data are obtained from the solution to the forward problem (expressions [[H) 
with known values of the loop parameters, so that the numerical inversions are performed under 
controlled conditions. In addition, this approach enables us to test the performance of the three 
different prior probability distributions for the density contrast discussed in Section [X2l 

We have considered a model coronal loop with the following physical parameters: ^ = 5, 
IIR = 0.25, and tm = 150 km s '. Forward modeling, according to the theory of resonantly 
damped kink waves under the thin tube and thin boundary approximations, predicts a period P = 
232.4 s and a damping ratio Td/P = 3.8. The period and the damping ratio are then used as data for 
the inversion and the ability of the code to recover the model coronal loop properties is analyzed, 
for different prior distributions and varying ranges for their definition. In all our simulations, 
variances of crp = O.IP and cr^^, = O.lTd are considered. 

Figures[2]and[3]show the inversion results obtained with these data, for fixed prior distributions 
for the transverse inhomogeneity length scale and internal Alfven travel time. For the density 
contrast, the three different prior types defined in Section IX^ are used in the range ^ e [1.5, 20]. 
In all three cases, the resulting Markov chains closely follow the same regions outlined by the 
analytic inversion curve. The distribution and density of the elements in the chain - which directly 
determine the properties of the posterior - differ in different regions of the parameter space. For 
a uniform prior distribution in density contrast, the elements in the chain are roughly uniformly 
distributed along the inversion curve (Figure [2^). The relevant quantitative information from the 
inversion comes from the analysis of the posterior probability distribution functions. Since the 
elements of the Markov chains are samples from the full posterior it is easy to divide the range 
for a given parameter in a series of bins and count the number of samples falling within each 
bin. By computing the marginal posterior for each parameter, we see that the density contrast 
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cannot be constrained (see Figure [3^). However, the projection of the chain onto the {l/R, tai)- 
plane indicates that a well defined posterior probability distribution can be obtained for these two 
parameters. Figures [Sj? and c show that, indeed, a very good estimation of these two parameters is 
possible. When a Jeffreys type prior for the density contrast is used, the elements in the chain are 
not uniformly distributed. Their density appears to be larger at positions in the parameter space 
that are close to the model coronal loop properties (Figure |2b)- The marginal posterior for the 
density contrast, however, does not show a well defined probability distribution (Figure |3^). The 
remaining two parameters, though, are well constrained (Figures |3]d and c). Finally, we introduce 
a hypothetical density contrast measurement for our model coronal loop and use a Gaussian prior 
distribution for the density contrast centered around that measurement. Figure|2fc shows that in that 
case, the inversion produces a Markov chain whose elements are closely packed around the correct 
loop properties. The marginal posteriors in this case (Figure |3]) show well-defined Gaussian-like 
distributions for the three parameters of interest. 

These results indicate that, in general, the Bayesian inversion scheme will enable us to con- 
strain two out of the three parameters of interest, the transverse inhomogeneity length-scale and 
the internal Alfven travel time. The density contrast, in general, cannot be constrained. If ad- 
ditional information about the internal and external densities is available from observations, this 
information can then be used to fully constrain the three unknowns. From the obtained posterior 
distributions in Figure [31 the median and the computation of the 68% confidence level can be used 
to obtain estimates with error bars for the inferred parameters. 

We found that the obtained results are largely insensitive to the variation of the range as- 
sumed for the internal Alfven travel time, as long as this range is sufficiently wide so it encloses 
the tails of the resulting posterior distribution. As for the density contrast, additional inversions 
were performed for different values of the upper limit in the prior distribution for this parameter, 
considering the three prior types. The results indicate that some prior types perform better than 
others. Figure S] displays the joint probability distribution for the internal Alfven travel time and 
the transverse inhomogeneity in the form of contours for the 68% and 95% confidence levels, for 
some illustrative cases. In all plots the known coronal loop parameters are indicated by a symbol. 
Numerical values for all performed inversions are given in Table [T] 

When a uniform prior for density contrast is considered, increasing the upper limit of the 
considered range (Figures |4^-c) produces a displacement of the joint marginal posterior, in such a 
way that the inversion underestimates the transverse inhomogeneity length scale and overestimates 
the Alfven travel time. Although the real values are for most of the cases enclosed within the 
error bars, we can see that for example, for the largest considered range (Figure |4t), the real value 
is already outside the 68% confidence level given by the inversion. Another effect that can be 
appreciated looking at Figures S^-c is the shrinking of the joint posterior for IjR, that results in 
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smaller error bars for this parameter (see also Table[T])- This result cannot be taken as an improved 
inversion, as is solely due to the fact that by extending the range in ^, the uniform prior gives more 
weight to elements in the chain with large density contrast, hence lower transverse inhomogeneity 
length scales and larger Alfven travel times are obtained. As a matter of fact, the integral of the 
prior distribution in Equation ([7]) in the range ^ 6 [10, 100] is about ten times the same integral in 
the range ^ e [1,10]. The optimal results for the inversion using a uniform prior in density contrast 
are obtained when the upper limit in ^ is around 2 or 3 times the real density contrast value. The 
problem is that in a real application we can hardly know the real density contrast value. 

A better performance and independence of the results on the prior information is obtained 
when considering a Jeffreys type prior for density contrast. Since this distribution assigns de- 
creasing probability for increasing contrast, the integrals of the prior in Equation ^ in the ranges 
^ G [1, 10] and ^ e [10, 100] are of the same order. Figures Uli-f show the results for this prior. 
It can be appreciated that the obtained marginal posteriors are almost independent of the assumed 
range on the prior information, once the upper limit for ^ is sufficiently large, and that optimal re- 
sults that recover the known coronal loop parameters are obtained (Table [I])- For this reason, even 
if the density contrast itself cannot be constrained, the use of a Jeffreys prior in density contrast is 
found to be appropriate to perform the inversion of the remaining two parameters in real coronal 
loops, in the case information about their density is lacking or uncertain. 

When density measurements are available, the Gaussian prior defined in Equation (ITOl) can 
be used. We have tested the performance of the inversion using this type of prior. As Figures [3l 
H^-i, and Table [T] show, in this case a full determination of the three unknowns is obtained, that 
remarkably well recovers the known coronal loop parameters. In addition, the inversion results are 
independent on the underlying assumptions of the a priori ranges in density contrast. Caution is 
called to the fact that, in this case, the reliability of the inferred parameters is entirely dependent 
on the reliability of the density contrast measurement. 

Additional inversions under controlled conditions were performed for different combinations 
for the model coronal loop properties with, e.g., ^ = 2.5, 10, 15 and l/R = 0.5,0.85. The obtained 
results do not change in a qualitative way those described above. 



4.2. Coronal loop oscillations 

Once the performance of our inversion method is evaluated, using synthetic dat a, we have ap 



plied t h e Bayesian inver s ion sc h eme to the 1 1 loop o s cillati o n events analyzed in e .g .JOfman & Aschwanden 



(2002); iGoossens et all (120021) : lAschwanden et all (l2003h : lArregui et all (120071) : iGoossens et al 



(|2008|) . Table |2] displays the main oscillation properties for these events and the results from the 
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application of analytic and Bayesian inversion techniques. Estimated values correspond to the 
median of the obtained distribution and errors are given at 68% confidence level. 

Observed period and damping ratios are used as data in our inversion code. Measurement 
errors in o bserved periods and darn ping times for oscillating coronal loops are lacking in the 
literature. lAschwanden et al.l (|2002h presented the most complete catalog of events with mea- 
surements of periods and damping times for 26 oscillating coronal loops. Uncertainties up to 
40 -60% can be found in those ev ents for the measured periods and damping times. The analysis 
by IVan Doorsselaere et al.l (l2007h drastically reduces the uncertainties to about 1 % in period and 
3% in damping time. It is not clear the meaning of such small errors on periods and damping 
times, with uncertainties that are shorter than the exposure time. In our analysis, we consider 
Gaussian errors of 10% on both quantities. Uniform prior distributions are considered for the 
transverse inhomogeneity length-scale, in the range l/R e [0, 2], and for the internal Alfven travel 
time, in a range determined by the oscillation period. We found that the ranges tai 6 [1, 400] and 
TAi 6 [1, 800] easily accommodate the resulting posterior distributions for the shorter and longer 
period oscillation events in Table [2l respectively. As for the prior information in density contrast, 
based on our assessment in Section WA\ we have performed the inver sions using both a J effrey s 
prior and a Gaussian prior, centered on the measurements reported by lAschwanden et al.l (|2003[) . 
for the very same 1 1 loop oscillation events. For the Jeffreys prior the minimum, ^min, is computed 
using the interval for ^ in expressions Q. The maximum value is set to ^max = 50, based on our 
assessment with synthetic results (Table [I])- For the Gaussian prior, besides t he estimated values 



we al so use the computed error bars (see column "Density ratio" in Table 3 by lAschwanden et al 



20031) . to implem ent the standard deviation in Equation (fTUl ) for each case. Uncertainties in density 
contrast listed by lAschwanden et al.l (|2003|) can be as high as 50%, because of the indirect methods 
followed to obtain those estimates. As the inversions using Gaussian prior information in density 
contrast are independent of the assumed ranges, we set the same range as before, ^ 6 [^min. 50]. The 
use of a uniform prior in density contrast is discarded, because of its relatively poor performance 
in front of the Jeffreys prior, in the absence of observational information on this parameter. 

Table |2] first lists the analytically obtained intervals for internal Alfven travel time computed 
using the corresponding expression for y = T/^i/P in expressions Q. This is all the information 
we can get from that inversion scheme, unless some additional assumption is made. Estimates for 
the Alfven travel time using both Bayesian inversions, with Jeffreys and Gaussian priors in density 
contrast, very well accommodate within the analyticall y obtained intervals. These es timates are di- 
rectly linked to the oscillation period, but in contrast to lNakariakov & Of man! (|200l|) . the inversion 
now makes use of the information contained in both the period and the damping time. In addition, 
the Bayesian inversion enables us to compute error bars for this parameter, which are consistently 
calculated from the uncertainties on data. Error bars in the determination of the Alfven travel time 
are found to be rather symmetric. 
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We next consider the determination of the transverse inhomogeneity length-scale. Regard- 
less of the prior type used in the inversion, our estimates for l/R are closely related to the observed 
damping rate, P/t^. Note that the larger this quantity, the stronger the damping by resonant absorp- 
tion is and the larger the obtained inhomogeneity length-scales are. Error bars in the determination 
of l/R are found to be rather asymmetric, with larger upper errors. The reason is that the posterior 
for l/R has a more elongated tail towards the right of the maximum, because of the contribution of 
samples with low values of the density contrast (see the projection onto the [^, //7?] -plane in Fig- 
ure[2b)- The transverse inhomogeneity length-scale is a quantity that is meant to capture the radial 
variation of the local Alfven frequency, which is unknown. Our modeling assumes a sinusoidal 
density variation. This being an assumption, it is not surprising that estimates for l/R are subject 
to uncertainties. By comparing the estimates of IjR obtained with Jeffreys and Gaussian priors, we 
find that the numerical estimate considerably differs in some cases, e.g., for loops #4, #5, #9, #10, 
#1 1 . When the errors bars are taken into account they are rather compatible. More importantly, the 
concordance between damping ratio and inferred transverse inhomogeneity length- scale holds in 
both cases, in such a way that an ordering of inhomogeneities can be established by just following 
an ordering of damping ratios. 



For the analyzed loop oscillation events. I Aschwanden et al.l (|2003b additionally provide esti 



mates for the transverse inhomogeneity length-scale, since the radius {R) and the skin depth (/) are 
measured (see their Table 2). These measurements lead to values of IjR in between 0.75 and 0.96, 
which in view of the values displayed on our Table [2l would point to rath er more transversely inho - 



mogeneous loops, than those obtained from our inversions. Notice that I Aschwanden et al.l (|2003|) 



have typical uncertainties of 20% (up to 50% in some cases) on their estimates of /. 

Density contrasts in the las t column basically recove r the input mean and standard deviation 



taken from the data measured by lAschwanden et al.l (l2003b and incorporated in the Gaussian prior. 



Their value is in the fact that they enable us to collapse the inversion chain, as explained in Sec- 
tion UTTl and more accurately determine the remaining two parameters, thus producing the smaller 
error bars in IjR, in comparison with the Jeffreys prior results. As mentioned above, the reliability 
of the Gaussian inversion, especially in the determination of the transverse inhomogeneity length- 
scale, is entirely dependent on the reliability of the density contrast measurement. For this reason, 
it is reassuring that estimates of IjR using both prior types - the less informative Jeffreys prior and 
the more informative Gaussian prior - result in similar values for most of the cases. On the other 
hand, because density measurements in the solar corona are challenging, we might consider that 
the inversion results using the less informative Jeffreys prior are more reliable that the ones that 
make use of uncertain information on density contrast. If that is the case, the differences in the 
determination of l/R when comparing both Bayesian inversions could be ascribed to inaccurate 
density contrast measurement. 
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5. Discussion and conclusions 



We have applied a Bayesian parameter inversion technique to the determination of unknown 
physical parameters in transversely oscillating coronal loops. The model considers coronal loops 
as one-dimensional cylindrically symmetric density enhancements with a radial density variation. 
The forward problem reduces to an algebraic set of equations for two observables as a function 
of three parameters. A Markov Chain Monte Carlo algorithm is used to sample the posterior 
probability distribution for the density contrast, the transverse inhomogeneity length-scale, and the 
internal Alfven travel time. We find that the latter two parameters can be very well determined, 
while the density contrast cannot in general be constrained. 

Density contrast estimations are challenging in both observations and MHD seismology in- 
versions. By using synthetic data to perform the inversion under controlled conditions, different 
prior information for the density contrast has been tested. A uniform prior distribution is unable 
to correctly infer the parameters of interest, unless some a priori information on density contrast is 
known, since it produces biased results when extended ranges are considered. On the other hand, 
a Jeffreys type prior is more appropriate to perform the inversion when information on densities is 
lacking or uncertain. If additional information on density contrast is available from observations, 
the three unknowns can be very well determined, by using a Gaussian prior centered on the ob- 
served density contrast value. The reliability of the inversion then depends on the reliability of the 
density contrast measurement. 



Th e advantage of the propose d technique with respect to Nakariakov & Ofmanl(|200l|) : lGoossens et al 



(120021) : lAschwanden et al.l (l2003h is that it uses the information contained in both the period and 
the damping rate of oscillating coronal loops. In addition it makes no assumption on the particular 
value of th e unknown paramete r s. The advantage with respect to the analytic and numerical in- 
versions by lArregui et al.l (|2007h : iGoossens et al.l (|2008|) is that, first, it enables us to constrain the 
transverse inhomogeneity length-scale and the density contrast. Second, the method incorporates 
confidence levels and error bars consistently calculated from the uncertainties on observed wave 
properties. Bayesian inference gives a measure of degree of belief about a proposition (set of pa- 
rameters) given some conditionals (observed data), once a probability model is set up. In that sense 
it is appropriate to perform inversio ns by making the few er p ossible assumptions in the compar- 
ison to observations, as proposed by lArregui et al.l (120071) and iGoossens et al.l (|2008b . Our study 
is based on a theoretical model that assumes the classic straight magnetic flux tube representation 
for coronal loops and resonant absorption as damping mechanism. If any of these assumptions 
deviates from reality, then this will int roduce errors into the inversion presented in this paper. For 
instance, iDe Moortel & Pascod (|2009|) have shown that the magnetic field strength inside a three- 
dimensional numerical coronal loop model may substantially differ from that inferred using classic 
seismology inversions in straight magnetic flux tubes. 
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The presented inversion technique can be applied to the determination of physical parameters 
in other magnetic and plasma structures of the solar atmosphere that display oscillatory dynam- 
ics. One example are transversely oscillating pror ninence threads, p rovided a large part of the 



magnetic flux tube is filled with dense plasma (see ISoler et alj|2010t) . In this case, and because 
of the large density contrast typical of prominence plasmas, with internal densities two orders of 
mag nitude larger than coronal densities, this parameter becomes irrelevant to the inversion pro- 



cess (lArregui & Ballesteij|2010l) . If we apply our algorithm to a case with P = 3 min and t^ = 9 
min (with 10% error), by considering large density contrast limits in the prior information for this 
parameter, we find estimates of l/R =0.22;^° J", tai =129;|;[3 j, in excellent agreement with the ana- 
lytic inversion approximation of l/R = 0.21, tai = 126.9. Although it may seem that in this case, 
as in the coronal loop inversion with Gaussian prior information, a two parameter-two observable 
problem seems to be fully constrained, we should be aware that we are dealing with incomplete 
information, because of uncertainties from measured quantities. The advantage of our method then 
lies in the consistent propagation of errors. 

One of the most appealing applications of Bayes' rule is model comparison, where the evi- 
dence in Equation ([3]) plays the fundamental role. Provided difli'erent damping mechanisms could 
be properly parameterized with similar physical model parameters and observables, a Bayesian 
model comparison could be performed in order to shed light on the still highly debated mechanism 
that produces the damping of transverse oscillations in coronal structures. 

LA. acknowledges the funding received under the project AYA2006-07637 by Spanish Min- 
isterio de Ciencia e Innovacion (MICINN) and FEDER Funds. A.A.R. acknowledges the financial 
support of the Spanish Ministerio de Ciencia e Innovacion (MICINN) through project AYA2010- 
18029 (Solar Magnetism and Astrophysical Spectropolarimetry) and CONSOLIDER INGENIO 
CSD2009-00038 (Molecular Astrophysics: The Herschel and Alma Era). This work was moti- 
vated by discussions held at the II Reunion Espanola de Fisica Solar y Heliosferica in the frame- 
work of RIA (Red de Infraestructuras de Astronomia) and funded by the Direccion General de 
Planificacion y Coordinacion of the Spanish MICINN, under the action CAC-2007-48. We are 
grateful to Jaume Terradas and Ramon Oliver for valuable comments. I. A. dedicates this paper to 
the memory of Ezequiel Hierro Palacio. 



REFERENCES 

Adamakis, S., Walsh, R. W., & Morton- Jones, A. J. 2010, Sol. Phys., 262, 1 17 
Andries, J., Arregui, I., & Goossens, M. 2005, ApJ, 624, L57 



17 



Andries, J., van Doorsselaere, T., Roberts, B., Verth, G., Verwichte, E., & Erdelyi, R. 2009, 
Space Sci. Rev., 149, 3 

Arregui, L, Andries, J., Van Doorsselaere, T., Goossens, M., & Poedts, S. 2007, A&A, 463, 333 

Arregui, L, & Ballester, J. L. 2010, Space Sci. Rev., in press, doi: 10. 1007/sl 1214-010-9648-9 

Aschwanden, M. J., De Pontieu, B., Schrijver, C. J., & Title, A. M. 2002, Sol. Phys., 206, 99 

Aschwanden, M. J., Fletcher, L., Schrijver, C. J., & Alexander, D. 1999, ApJ, 520, 880 

Aschwanden, M. J., Nightingale, R. W., Andries, J., Goossens, M., & Van Doorsselaere, T. 2003, 
ApJ, 598, 1375 

Asensio Ramos, A. 2011, ApJ, 731, 27 

Asensio Ramos, A., Martinez Gonzalez, M. J., & Rubino-Martin, J. A. 2007, A&A, 476, 959 

Bayes, M., & Price, M. 1763, Royal Society of London Philosophical Transactions Series I, 53, 
370 

De Moortel, I., & Pascoe, D. J. 2009, ApJ, 699, L72 

Edwin, P M., & Roberts, B. 1983, Sol. Phys., 88, 179 

Goossens, M. 2008, in lAU Symposium, Vol. 247, lAU Symposium, ed. R. Erdelyi & 
C. A. Mendoza-Briceiio, 228-242 

Goossens, M., Andries, J., & Arregui, I. 2006, Royal Society of London Philosophical Transactions 
Series A, 364, 433 

Goossens, M., Andries, J., & Aschwanden, M. J. 2002, A&A, 394, L39 

Goossens, M., Arregui, I., Ballester, J. L., & Wang, T. J. 2008, A&A, 484, 851 

Goossens, M., Terradas, J., Andries, J., Arregui, I., & Ballester, J. L. 2009, A&A, 503, 213 

Gregory, P. C. 2005, Bayesian Logical Data Analysis for the Physical Sciences: A Comparative 
Approach with 'Mathematica' Support (Cambridge University Press) 

Ireland, J., Marsh, M. S., Kucera, T. A., & Young, C. A. 2010, Sol. Phys., 264, 403 

Marsh, M. S., Ireland, J., & Kucera, T. 2008, ApJ, 681, 672 

Metropolis, N., Rosenbluth, A. W., Rosenbluth, M. N., Teller, A. H., & Teller, E. 1953, 
J. Chem. Phys., 21, 1087 



18 



Nakariakov, V. M., & Ofman, L. 2001, A&A, 372, L53 

Nakariakov, V. M., Ofman, L., DeLuca, E. E., Roberts, B., & Davila, J. M. 1999, Science, 285, 
862 

Neal, R. M. 1993, Probabilistic Inference Using Markov Chain Montecarlo Methods, Tech. Rep. 
0560, Dept. of Statistics, University of Toronto 

Ofman, L., & Aschwanden, M. J. 2002, ApJ, 576, L153 

Roberts, B., Edwin, R M., & Benz, A. O. 1984, ApJ, 279, 857 

Ruderman, M. S., & Erdelyi, R. 2009, Space Sci. Rev., 149, 199 

Soler, R., Arregui, I., Oliver, R., & Ballester, J. L. 2010, ApJ, 722, 1778 

Tandberg-Hanssen, E. 1995, The nature of solar prominences (Dordrecht ; Boston : Kluwer, 
cl995.),358 

Trotta, R. 2008, Contemporary Physics, 49, 71 

Uchida, Y. 1970, PAS J, 22, 341 

Van Doorsselaere, T., Andries, J., Poedts, S., & Goossens, M. 2004, ApJ, 606, 1223 

Van Doorsselaere, T., Nakariakov, V. M., & Verwichte, E. 2007, A&A, 473, 959 

Verth, G., Erdelyi, R., & Jess, D. B. 2008, ApJ, 687, L45 

Verwichte, E., FouUon, C, & Nakariakov, V. M. 2006, A&A, 452, 615 

Zaqarashvili, T. V. 2003, A&A, 399, L15 



This preprint was prepared with the AAS lATj^ macros v5.2. 



19 











2 



4 6 



8 



10 



Fig. 1 . — Example of Uniform, Jeffreys, and Gaussian prior probability distributions for the density 
contrast in the range ^ e [1, 10]. For the Gaussian distribution, a mean value ^ = 5 with standard 
deviation cr^ = 0.5 has been considered. 
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Fig. 2. — Three-dimensional view (and projections onto the different planes) of the converged 
Markov chains in the (^, l/R, tai) parameter space for a synthetic coronal loop with ^ = 5, IjR = 
0.25, and Tai = 150 km s~^ The inversions were performed using uniform prior distributions for 
l/R 6 [0, 2] and tai 6 [1, 400] and three different prior distributions for ^ e [1.2, 20]: (a) uniform 
prior; (b) Jeffreys prior; and (c) Gaussian prior centered dX ^ = 5 and variance a^ = O.IZ;". In 
all plots, solid-lines correspond to the analytic inversion curve and variances of crp = O.IP and 
ar^ = 0. iTd have been used for the oscillation period and damping time. 
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Fig. 3. — One-dimensional marginalized posterior distributions for the density contrast (a), the 
transverse inhomogeneity length- scale (b), and the internal Alfven travel time (c), corresponding 
to the inversions displayed in Figure [21 The solid, dashed, and dotted lines correspond to the three 
different prior probability distributions for density contrast, as indicated in Figure [IJ 
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Fig. 4. — Joint two-dimensional posterior distributions for the transverse inhomogeneity length- 
scale and the internal Alfven travel time for different combinations of prior types and ranges of 
variation for the density contrast. From left to right, (a)-(c) uniform prior; (d)-(f) Jeffreys prior; 
and (g)-(i) Gaussian prior with H( = 5 and a^ = 0.5 . From top to bottom, (a)-(g) ^ 6 [1.2, 10]; 
(b)-(h) ^ e [1.2,20]; and (c)-(i) ^ e [1.2,50]. The outer boundaries of the light grey and dark 
grey shaded regions indicate the 95% and 68% confidence levels. Symbols indicate the synthetic 
coronal loop properties, l/R = 0.25, and tai = 150 km s"^ 
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Table 1 . Inversion of synthetic data for different priors types and ranges of variation in density 

contrast. 



Prior type 



Uniform 



Jeffreys 



Gaussian 



^-range 



IjR 



TAi (S) 



HR 



TAi (S) 



^ 



IjR 



TAi (S) 



1.2-10 
1.2-15 
1.2-20 
1.2-30 
1.2-40 
1.2-50 
1.2-60 
1.2-100 



36+"™ 

22+°" 
21+009 

19+0.06 
19+0.05 
19+""4 
19+""4 




0.52 
0.53 
0.51 
0.51 
0.51 
0.50 
0.49 
0.52 
0.52 
0.52 
0.51 
0.51 
0.49 
0.51 
0.49 
0.48 



26+""4 
26+""4 
25+""4 
25+°"4 
25+°"4 
26+°"4 
25+°"4 
26+""4 



152 4+15.4 

152 1+15° 

151 1+15-^ 

152 0+15-^ 
152.0!}5-3 

151.9!}^-^ 
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Table 2. Analytic and Bayesian inversion results for the analyzed loop oscillation events. 



Oscillation properties 



Analytic 



Inversion results 



Bayesian Jeffreys 



Bayesian Gaussian 



# P(S) Td(s) P/Td 



TAi (S) 



TAi (S) 



HR 



TAi (S) 



HR 



1 


261 


870 


0.30 


145-177 


161 5+22 2 


36+"-2^ 


169.4:;^:^ 


2 


265 


300 


0.88 


163-182 


169.9!20-^ 


92+"-47 


167.1!!^:^ 


3 


316 


500 


0.63 


189-217 


199.4:250 


76+"'^' 

^•'"-0.28 


196.8:j^-4 


4 


277 


400 


0.69 


168-189 


iie.i-j^] 


73+"-^3 


167.2_+I^:^ 


5 


272 


849 


0.32 


151-187 


173.2+21-1 


34+0-26 

^••^ -0.11 


i59.7_+;^:^ 


6 


522 


1200 


0.44 


304-359 


329 7+43-8 


49+0-39 


319.9_+1^:^ 


7 


435 


600 


0.73 


267-299 


281.3!33-i 


74+"-4i 

'-'•'^-0.20 


290.9!j^-4 


8 


143 


200 


0.72 


90-98 


90.9!}2-o 


76+°-53 

'-'•'"-0.23 


93.8:1^:^ 


9 


423 


800 


0.53 


247-291 


265 6+3^-2 


64+"'^^ 


290.5:;^;4 


10 


185 


200 


0.93 


117-126 


1 19 2+i'*-^ 


94+0-48 


1 14 4+17-4 


11 


390 


400 


0.98 


245-270 


250.5^29-6 


A gQ+0.54 
^•^^-0.28 


221 5+"-'* 



30+""^ 

"-'"-0.04 

1 01+°-''' 

77+0-43 

"•' -0.19 

1 05+°'^3 

i.u^_Q28 

0.58t- 
0.59tf3 
64+°" 

"•"^-0.09 

69+0-11 
"•"'^-o.io 

41+0.07 

1 9 1 +0.24 
^•^^-0.20 

1 69+°-'^ 



4 QQ+0.50 

3 76+0-64 

o ro + 1.88 

2 56+"-''^ 

1 O+0.75 

9 Q7+0.94 
^•^^'-0.91 

6.98!1- 

r rr+0.94 

J.JJ Qgg 

1 -i 4+3.50 
-^' -3.80 

3.08-- 
2. 10-- 



